Note the change in the urban center. The general pattern of blackouts are visible in the form of orange encroaching into the yellow area. This is apparent despite the overall increase in light intensity across the area.
Part 2: Identifying impacted homes and census tracts based on changes in night light intensity.
Homes impacted by blackout
Code
# use h_light_0207 to make a Houston area vector for later usehouston <-st_as_stars(h_light_0207)# vectorizehouston <-st_as_sf(houston)# check / convert Houston cookie cutter data CRSif(crs(houston) =="EPSG: 3083"){print("CRS correct")} else{warning("Updated CRS to EPSG 3083") houston <-st_transform(houston, "EPSG: 3083")}# dissolve the grid into cookie cutterhouston <-st_union(houston)# Calculate light intensity differencelight_diff <- light_02_07 - light_02_16summary(light_diff)# Calculate light intensity difference## define reclass matrixrcl <-matrix(c(-20900, -200, 1,# group 1 ranges from -20900 (min) to -200-201, 3775, NA), # group 2 ranges from -200 to 3776 (max)ncol =3, byrow =TRUE)# use reclass matrix to reclassify light_diff rasterblackouts <-classify(light_diff, rcl = rcl)# change reclasssed values into factorsvalues(blackouts) <-as.factor(values(blackouts))summary(blackouts)# Crop blackout raster to Houston areablackouts_cropped <- terra::crop(blackouts, bbox)summary(blackouts_cropped)# Vectorize light raster## first, convert it to a stars objectblackouts_stars <-st_as_stars(blackouts_cropped)# have a look at itsummary(blackouts_stars)# vectorize!!!blackouts_vector <-st_as_sf(blackouts_stars)# check / convert blackout vector CRSif(crs(blackouts_vector) =="EPSG: 3083"){print("CRS correct")} else{warning("Updated CRS to EPSG 3083") blackouts_vector <-st_transform(blackouts_vector, "EPSG: 3083")}
Code
# load county data as a basemapcounties <-st_read(here("data", "tl_2023_us_county", "tl_2023_us_county.shp"))# check / convert county vector CRSif(crs(counties) =="EPSG: 3083"){print("CRS correct")} else{warning("Updated CRS to EPSG 3083") counties <-st_transform(counties, "EPSG: 3083")}# crop the county data to the Houston areahouston_counties <-st_filter(x = counties, y = houston,.predicate = st_intersects)# Make a 200 m buffer around roads## check for meters as unitsst_crs(osm_roads)$units # buffer the roads by 200mhighway_buffer <-st_buffer(osm_roads, dist =200)# combine the buffer polygonshighway_buffer_union <-st_union(highway_buffer)## Exclude highway buffer from blackout maskblackouts_mask_no_highways <-st_difference(blackouts_vector, highway_buffer_union)# Spatial geometry filtering to count homes impacted by blackoutsblackout_homes <-st_filter(x = osm_buildings, y = blackouts_mask_no_highways, .predicate = st_intersects)
## Crop census data to the Houston areacensus_income <-st_filter(x = census_income, y = houston,.predicate = st_intersects)# filter down to the affected tractsblackout_tracts <-st_filter(x = census_income, y = blackouts_mask_no_highways,.predicate = st_intersects)# filter down to the unaffected tracts ## (st_difference isn't working so this is my work-around)# find the indices of blackout polygonsblackout_indices <-st_intersects(census_income, blackouts_mask_no_highways)# filter by those indicesunaffected_tracts <- census_income[!sapply(blackout_indices, length) >0, ]
Homes likely to have experiencced the blackout are depicted as pink polygons in the above figure. An estimated 81914 homes experienced blackouts.
Here the census tracts likely to have experienced the blackout are visible as pink polygons overlaid on the Houston area.
Part 3: Comparing the median household income among impacted and un-impacted census tracts.
Code
# Add a grouping variable to your data beforehand if neededblackout_tracts$status <-"Impacted"unaffected_tracts$status <-"Unaffected"# Combine the datacombined_tracts <-rbind(blackout_tracts, unaffected_tracts)plot_blackout_tracts <-ggplot(data = combined_tracts, aes(x = B19013e1, fill = status)) +geom_histogram(data = blackout_tracts, alpha =0.5, bins =30) +geom_histogram(data = unaffected_tracts, alpha =0.5, bins =30) +scale_fill_manual(name ="Census Tract Status", # Legend titlevalues =c("Impacted"="#fd007b", "Unaffected"="#3b00ab"),labels =c("Impacted by Blackout", "Not Impacted")) +# Legend labelslabs(x ="Median Household Income ($)",y ="",title ="Distribution of Median Incomes for Census Tracts Impacted by Blackouts") +theme_bw() +theme(panel.background =element_rect(fill ='transparent'),plot.background =element_rect(fill ='transparent', color =NA))ggsave(here("outputs", "plot_blackout_tracts.png"), plot_blackout_tracts,height =45, width =7, dpi =400, bg ="transparent")
Figure Interpretation
Above are depicted wide spread of incomes impacted by the blackouts. Note that lower median household incomes are more highly represented than un-impacted tracts.
Discussion
In the above analysis and accompanying figures, the impacts of the storm are made apparent as a visible change in night light intensity in the hoston metro area. An estimated 81914 homes experienced blackouts according to this analysis. Further, these results suggest that blackouts occured independently of income level, as indicated by the median household incomes across census blocks. This is somewhat surprising as one might expect households with higher income to overcome barriers and power their homes (and lights). This trend implies the pervasive nature of the blackouts without readily available alternate power sources.
Acknowledgements
This analysis and workflow was originally created by Dr. Ruth Oliver of the Bren School of Environmental Science & Management for the Masters of Environmental Data Science course, Geospatial Analysis and Remote Sensing. My initial work through of this analysis was conducted as part of a homework assignment for this class, and I later polished up this repository.
Data
The data used in this analysis is too large to be hosted on GitHub. Instead, download the data here as a zipped folder, unzip, and move into the R project manually.
Citations
Data
Link
Night Lights
https://ladsweb.modaps.eosdis.nasa.gov/
Open Street Map Roads
https://planet.openstreetmap.org/
Open Street Map Buildings
https://planet.openstreetmap.org/
Socioeconomic Data
https://www.census.gov/programs-surveys/acs
Citation
BibTeX citation:
@online{mitchell2024,
author = {Mitchell, Steven},
title = {Houston {Storm} {Blackouts}},
date = {2024-11-09},
url = {https://steven-mitchell.github.io/work-samples/houston-storm-blackouts},
langid = {en}
}